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Bose-Einstein condensates in standing waves: The cubic nonlinear Schrodinger 

equation with a periodic potential 
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We present a new family of stationary solutions to the cubic nonlinear Schrodinger equation with 
a Jacobian elliptic function potential. In the limit of a sinusoidal potential our solutions model a 
dilute gas Bose-Einstein condensate trapped in a standing light wave. Provided the ratio of the 
height of the variations of the condensate to its DC offset is small enough, both trivial phase and 
nontrivial phase solutions are shown to be stable. Numerical simulations suggest such stationary 
states are experimentally observable. 



The dilute-gas Bose-Einstein condesate (BEC) in the 
quasi-one-dimensional regime is modeled by the cubic 
nonlinear Schrodinger equation (NLS) with a poten- 
tial The various traps which are used to con- 
tain the BEC have spurred the solution of the NLS with 
new potentials [Qj5). BECs trapped in a standing light 
wave have been used to study phase coherence || and 
matter-wave diffraction [Q and have been predicted to 
have applications in quantum logic H and matter-wave 
transport (^] . Exact solutions have been obtained for the 
Kronig-Penney potential (!(]] and some researchers have 
used a Bloch function description |0j]. In this letter, we 
study new explicit solutions of the NLS with a Jacobian 
elliptic function potential. 

We consider the mean-field model of a quasi-one- 
dimensional repulsive BEC trapped in an external po- 
tential which is given by the nonlinear Schrodinger equa- 
tion | 



(1) 



In experiments, the trapping potential is generated by a 
standing light wave |J. As a model for such a potential 
we use the periodic potential 



V(x) = -V sn 2 (x,/c) 



(2) 



where sn(x, k) denotes the Jacobian elliptic sine func- 
tion |l2| with elliptic modulus < k < 1. In the limit 
k — > l - , V(x) becomes an array of well-separated hyper- 
bolic secant potential barriers or wells, while in the limit 

— > + it becomes purely sinusoidal. We note that for 
intermediate values (e.g. k = 1/2) the potential closely 
resembles the sinusoidal behavior and thus provides a 
good approximation to the standing wave potential gen- 
erated experimentally ||. 

We present stationary solutions in closed form and 
study their stability analytically and numerically. We 
begin by constructing solutions to Eq. ([!]) which have 
the form 



ip(x, t) = r(x) exp(z(0(x) — cot)) 



where 
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where B determines a mean amplitude and acts as a DC 
offset for the number of condensed atoms. The strength 
of the nonlinearity, which for the BEC is a function of 
both the atomic coupling and the number of condensed 
atoms, is determined by the parameters Vo + k 2 and B, 
as is apparent in the amplitude of the solutions given by 
Eq. (||). Note that if x is scaled so that V(x) undergoes 
only a single oscillation on the ring (in the limit k — » 1) 
the Jacobian elliptic potential provides a model of a sin- 
gle barrier or well (l3). For simplicity we focus on two 
special cases: (1) k arbitrary and trivial phase (c = 0), 
and (2) k = with nontrivial phase (c ^ 0) so that the 
solutions are trigonometric functions. 

Trivial Phase Case - In the limit of c = 0, the solutions 
given in Eqs. (||)-(fl]) reduce to 



sn(x, fc)expH(l + k 2 )t/2] (5) 

valid for Vb > — k 2 , and 

4>(x, t) = ^/-(V + k 2 ) cn(x, k) exp[-i(V + k 2 - ±)t], (6a) 

tftx, t) = ^/-(l + p) dn(x, k) exp[i(l + ^ - y )t] (6b) 

valid for Vb < — k 2 where cn(x, k) and dn(x, k) are Jaco- 
bian elliptic functions ll2| . These solution branches 
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FIG. 1. Regions of validity for trivial phase Jacobian 
elliptic solutions. The sn(x,k) and cn(x,k) branches are 
found to be unstable while the dn(x,k) branch is stable. 



are illustrated in Fig. 1 along with their stability prop- 
erties found from analytic calculations and observed in 
numerical simulation. 

We can prove that the dn(x, k) branch of solutions is 
linearly stable. To do so, we linearize about the dn(x, k) 
solution branch given by Eq. (pp) so that 



ip(x,t) = (4>o(x) + t)) exp(— tout) 



(7) 



where (j)o(x) exp(— iuot) is the exact solution given by 
Eq. (^|b) and (f>(x, t) <C 1 is a small perturbation to the ex- 
act solution. This leads to the following linearized eigen- 
value problems 



— A 2 R 
L+L-I = X 2 I 



(8a) 
(8b) 



where (j)(x,t) = (R + il) exp(i\t) is decomposed into its 
real and imaginary parts. The operators L_ and L + are 
both self-adjoint and periodic differential operators: 
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^(x) 2 - V{x)-uj (9a) 
3^ (^) 2 - V(x) -co. (9b) 



Thus the eigenvalue problem closely resembles that of 
soliton solutions of the NLS with the additional dif- 
ficulty of V(x) being a periodic potential. 

We are able to prove stability for the whole dn(x, k) 
branch of solutions and instability for the cn(x, k) and 
sn(x, k) branches near their emanation point Vb = — k 2 
(see Fig. 1). We note that stability or instability is proven 
by showing that the eigenvalue A 2 in Eq. (||) is positive 
or negative respectively. For the dn(x, k) branch of so- 
lutions, we find that L + is a positive definite operator. 
This allows us to construct its inverse. We then consider 
the related problem L_J = A 2 !^ 1 / which is self- adjoint 
with respect to the weighting operator L^ 1 . This then 
allows us to show that A 2 > so that the eigenvalues 
remain real, establishing stability. 




FIG. 2. Evolution of initially perturbed Jacobian el- 
liptic dn(x,/c) and sn(x,k) solutions. Panels (a)-(b) corre- 
spond to the stable dn(x,k) branch of solutions of Fig. 1 
with k = 1/2 and Vb — — 1. The wavenumber spectrum 
(77) remains constant and consists of the zero mode and 
wavenumber one. These reflect the DC offset and strength 
of oscillation, respectively. Panels (c)-(d) demonstrate the 
instability of the sn(x,fc) solution with k = and Vb = 1. 
Here the instability is seen to develop around wavenum- 
ber unity as predicted analytically. The cn(x,k) branch of 
solutions exhibits the same instability. 



In contrast to the dn(x, k) solutions, numerical experi- 
ments suggest that all cn(x, k) and sn(x, k) solutions are 
linearly unstable. Although we have been unable to show 
this for the whole solution branch, we have shown this 
perturbatively for Vb + k 2 « 1 using the same technique 
as outlined above for the dn(x, k) branch of solutions. In 
this limit, the dispersion relation for small disturbances 



near Vq + k 2 « 1 is 



(v 2 
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(10) 



where rj is the wavenumber of the disturbance and C is 
a non-zero constant. From this it follows that there is a 
band of unstable modes which occurs near 77 = 1. Numer- 
ical experiments show that this unstable band persists for 
all Vb and grows as Vb + k 2 becomes large. We empha- 
size that the stable dn(x, k) branch of solutions consists 
of oscillations about a non-zero mean value, while the 
unstable cn(x, k) and sn(x, k) branches have zero mean, 
suggesting that the offset has an important effect on the 
stability properties. Rigorous proofs of the linear stabil- 
ity calculations will be considered elsewhere [ p~3|| . 

Figure ^ illustrates numerical solutions of Eq. (|l|) for 
initial conditions consisting of the exact sn(x, k) and 
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where 



FIG. 3. Evolution of perturbed trigonometric solutions 
(k = 0) with nontrivial phase. In panels (a)-(b), B = 1 so 
that the solution is dn(x, /c)-like, leading to a stable evo- 
lution. The initial amplitude and phase are depicted in 
the top right. For B = 1/2, the solution is unstable, as 
seen in (c). Its initial phase and amplitude are depicted in 
panel (d). This indicates that a sufficiently high DC offset 
is required to stabilize the evolution. 



dn(x, k) solutions perturbed with a small amount of ini- 
tial stochastic white noise, along with the corresponding 
Fourier spectra. Figures ||(c) and ||(d) depict the evo- 
lution of a sn(x, k) initial condition perturbed by noise. 
For this simulation, Vq = 1 and k = so that the so- 
lution is a simple sinusoid (see Fig. 1). The solution is 
unstable, and diverges rapidly from the exact solution 
given in Eq. (||). From the plot of the Fourier spectrum, 
it is clear that the instability begins in a neighborhood 
of wavenumber rj = 1 as predicted from the linear stabil- 
ity analysis. Similar behavior is observed for the cn(x, k) 
branch of solutions. Figures §(a) and ||(b) depict the 
evolution of the stable dn(x, k) branch of solutions which 
is initially perturbed by white noise. In this simulation 
Vo = — 1 and k = 1/2. This solution branch is stable, and 
stays close to the exact solution for all time. Note that 
the Fourier spectrum is dominated by the zero mode, 
which determines the DC offset, and the wavenumber 
one, which determines the oscillation strength. 

The Nontrivial Phase, Trigonometric Limit - For 
k = 0, sn(x, 0) = sin(x). The governing evolution Eq. (|l|) 
reduces to 
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ipxx + H 2 ^ - Vo sin 2 0)^ . 



(11) 



Note that by the trigonometric identity 2 sin 2 (x) = 1 — 
cos(2x) this potential is sinusoidal. The solutions reduce 
to 

i/j(x, t) = \J V sin 2 x + B exp[i(6(x) - (1/2 + B)t)\ (12) 



tan(6(x)) = ±y/l + Vb/£tan(x) (13) 



determines the nontrivial phase provided that B > — Vo 
for Vo < and B > for V > 0. 

The results of the previous section imply that the triv- 
ial phase solutions without offset are unstable. Oscilla- 
tions about some mean value is qualitatively similar to a 
dn(x, k) solution, which suggests that solutions with suf- 
ficiently large B might be stable. Oscillations with \B\ 
small are qualitatively like cn(x, k) and sn(x, k) and are 
expected to be unstable. Numerical experiments confirm 
this. We note that a linear stability analysis in this case 
is more complicated than for the trivial phase case since 
the linearized operators do not decouple as in Eq. (||). 

In Fig. H we depict the evolution of a pair of initial con- 
ditions of the form given by Eq. (|l2|) plus initial white 
noise. Figures ||(c) and ||(d) depict the initial ampli- 
tude and phase along with the evolution of the density 
\ifj(x,t)\ 2 for the parameter values Vo = 1.0,5 = 0.5. It 
is clear from the graph that this stationary solution is 
unstable. Figures ||(a) and ||(b) depict the initial am- 
plitude and phase and the evolution of the density for 
parameter values Vo = 1.0, B = 1.0. In this case the so- 
lution appears stable. The solution with the perturbed 
initial condition stays close to the stationary solution for 
long times, which for our time scaling is longer than the 
lifetime of typical trapped BECs 

The numerical and analytical results imply that in 
order to obtain a stable condensate it is necessary to 
have solutions which are sufficiently in the nonlinear 
regime. To quantify this, we note that from Eqs. ( JTT| ) 
and ( |l2j) , the number of particles per well n is given by 
n = (J Q \ijj(x,t)\ 2 dx)/n; = Vo/2 + B. In the context of 
the BEC, and for a fixed atomic coupling strength, this 
means a sufficient number of condensed atoms per well 
n is required to provide a DC offset on the order of the 
potential strength. This ensures stabilization of the con- 
densate. 

To demonstrate the physical viability of our results, 
we perform numerical simulations consistent with recent 
experiments on the BEC ||. In particular, we consider 
a fixed number of condensed atoms with an initially con- 
stant DC profile. We adiabatically ramp up the periodic 
potential linearly from zero to a fixed potential strength 
Vq. In an experiment, the condensate initially has a triv- 
ial phase profile and the solution deforms as depicted in 
Fig. ||(a). This adiabatic process generates a stable solu- 
tion which is close to the dn(x, k) solution branch given 
by Eq. (||b) and depicted in Fig. ||(a). In contrast, an 
initial DC profile with a nontrivial phase profile is only 
stable provided the DC 
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FIG. 4. Stable growth of periodic condensate with DC 
offset of one with (a) initial trivial phase and (c) nontrivial 
phase. A sin 2 (x) potential was adiabatically grown from 
V = at t = to V = 0.5 at t = 50. Figure (b) shows the 
difference in the final state of the system given an initial 
trivial phase (black line) and initial nontrivial phase (gray 
line) . The development of phase structure on an initial lin- 
ear phase profile is depicted in (d). Here the initial linear 
phase profile has been subtracted out to give the initial 
(black line) and final (gray line) phase profiles. 

offset is sufficiently large in comparison to the potential 
strength Vq. For Vo = 0.5 and an initial DC offset of one, 
the adiabatically grown solution is stable with larger am- 
plitude fluctuations than those of the trivial phase case 
(Fig. ^b). As the condensate evolves, the initial linear 
phase profile is deformed as shown in Fig. [l](d), where the 
initial linear phase profile is subtracted out from the ini- 
tial and final phases. This phase deformation is necessary 
in order for the solution to remain stationary. In partic- 
ular, a linear phase profile induces a group- velocity shift 
in the direction of growing phase. This is in contrast to 
phase jumps or transition regions which can cause motion 
opposite the direction of a positive jump [||,||. The sta- 
tionary solutions generated here suggest that these two 
opposite actions effectively balance each other in order 
for the solutions to remain localized and stationary in 
their respective troughs. Finally, we note that for Vb = 1 
the potential is sufficiently strong so as to destabilize the 
adiabatically grown solution. The instability mechanism 
is similar to that observed in Fig. |](c). 

Conclusions - We have considered the repulsive NLS 
equation with a Jacobian elliptic function potential as 
a model for a trapped, quasi-one-dimensional Bose- 
Einstein condensate, and thereby produced a whole fam- 
ily of solutions. We have shown the stability of the 
dn(x, k) branch of solutions rigorously. In addition, we 
have used a perturbative argument to show the instabil- 
ity of the sn(x, k) and cn(x, k) solutions. In the sinusoidal 



limit we have provided numerical evidence that the non- 
trivial phase solutions are stable for sufficiently large off- 
set B. Stable trivial and nontrivial phase trigonometric 
solutions can be obtained through adiabatic ramping of 
the potential strength provided there is a sufficient DC 
offset in the number of condensed atoms. As this models 
a Bose-Einstein condensate trapped in a standing light 
wave, our results imply that sufficiently large number of 
condensed atoms (DC offset) are required to form a sta- 
ble, periodic condensate. 
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